Data Loading and Processing

First, we’ll load our dependencies and process the crash data:

# Source utility functions (which loads libraries)
source("R/utils.R")

# Read and process crash data
crashes <- load_crash_data("data/raw/crashes/traffic_crashes.csv")
crashes_sf <- convert_to_sf(crashes)

# Read and process fatality data
crash_fatalities <- load_fatality_data("data/raw/crashes/crash_fatalities.csv")
crash_fatalities_sf <- convert_to_sf(crash_fatalities, "Longitude", "Latitude")

# Combine crashes and fatalities
combined_crashes_sf <- bind_rows(crashes_sf, crash_fatalities_sf)

District and Street Processing

Next, we’ll process the senate district data and Lake Shore Drive information:

# Load senate district shapefile and filter for lakefront districts
senate_districts <- load_district_shapefile(
  "data/raw/shapefiles/senate/senate_districts.shp",
  new_name = "senate_district"
)
## Reading layer `senate_districts' from data source 
##   `/Users/mmclean/dev/r/ilga_district_crashes/data/raw/shapefiles/senate/senate_districts.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 59 features and 3 fields
## Geometry type: POLYGON
## Dimension:     XY
## Bounding box:  xmin: -91.51308 ymin: 36.9703 xmax: -87.4952 ymax: 42.50848
## Geodetic CRS:  NAD83
lakefront_districts <- senate_districts %>%
  filter(senate_district %in% c("6", "7", "13"))

district_6 <- senate_districts %>%
  filter(senate_district == "6")

# Spatial join crashes with district 6
crashes_with_district <- st_join(combined_crashes_sf, district_6, left = FALSE)

# Transform to meter-based CRS for accurate distance calculations
crashes_with_district <- st_transform(crashes_with_district, 3857)

# Query OSM for Lake Shore Drive
lake_shore_streets <- opq(bbox = CHICAGO_BBOX) %>%
  add_osm_feature(key = 'highway') %>%
  add_osm_feature(key = 'name', value = 'Lake Shore', value_exact = FALSE) %>%
  osmdata_sf() %>%
  pluck("osm_lines") %>%
  {rename(., street_name = name)} %>%
  st_transform(crs = 3857)

Crash Analysis

Process crashes within the Lake Shore Drive buffer:

# Create buffer around Lake Shore Drive and identify intersecting crashes
buffer <- st_buffer(lake_shore_streets, dist = BUFFER_DISTANCE_M)
intersects <- st_intersects(crashes_with_district, buffer, sparse = FALSE)

# Filter crashes within buffer and calculate counts
crashes_with_district$is_lakeshore_crash <- apply(intersects, 1, any)
lakeshore_crashes <- crashes_with_district %>%
  filter(is_lakeshore_crash) %>%
  mutate(
    crash_count = 1,
    fatality_count = ifelse(!is.na(FATALITY_VICTIM) & FATALITY_VICTIM != "", 1, 0)
  )

# Create fatality subset
lakeshore_fatalities <- lakeshore_crashes %>%
  filter(fatality_count > 0)

Crash Impact Visualizations

Fatalities Over Time

# Prepare fatalities data
fatalities_by_year <- lakeshore_fatalities %>%
  st_drop_geometry() %>%
  mutate(year = lubridate::year(CRASH_DATE)) %>%
  group_by(year) %>%
  summarize(total_fatalities = sum(fatality_count, na.rm = TRUE))

# Get min and max years for consistent x-axis
min_year <- min(fatalities_by_year$year)
max_year <- max(fatalities_by_year$year)

# Create the fatalities visualization
ggplot(fatalities_by_year, aes(x = year, y = total_fatalities)) +
  geom_col(fill = "#9370DB", width = 0.7) +
  geom_text(
    aes(label = total_fatalities),
    vjust = -0.5,
    size = 3.5,
    family = "Arial"
  ) +
  labs(
    x = "Year",
    y = NULL,
    title = "Traffic Fatalities Along Lake Shore Drive in District 6",
    subtitle = "Annual count of traffic-related fatalities",
    caption = "Data: Chicago City Data Portal \"Traffic Crashes - Vision Zero Chicago Traffic Fatalities\" Retrieved Feb. 6, 2025
Includes crashes reported within 100ft of DuSable Lake Shore Drive"
  ) +
  scale_x_continuous(
    breaks = seq(min_year, max_year, by = 1)
  ) +
  scale_y_continuous(
    expand = expansion(mult = c(0, 0.15)),
    breaks = function(x) seq(0, ceiling(max(x)), by = 1)
  ) +
  theme_classic() +
  theme(
    plot.title = element_text(face = "bold", size = 12, margin = margin(b = 10)),
    plot.subtitle = element_text(size = 10, color = "#666666"),
    plot.background = element_rect(fill = "#FFFCF7"),
    panel.background = element_rect(fill = "#FFFCF7"),
    panel.grid.major.y = element_line(color = "#f0f0f0", linewidth = 0.3),
    text = element_text(family = "Arial"),
    axis.line = element_line(color = "#666666"),
    axis.ticks = element_line(color = "#666666"),
    axis.title = element_text(color = "#666666"),
    axis.text = element_text(color = "#666666"),
    plot.caption = element_text(
      color = "#666666",
      size = 8,
      margin = margin(t = 10),
      hjust = 0
    ),
    plot.margin = margin(t = 15, r = 20, b = 15, l = 15, unit = "pt")
  )

Total Injuries Over Time

# Prepare injuries data
injuries_by_year <- lakeshore_crashes %>%
  st_drop_geometry() %>%
  mutate(
    year = lubridate::year(CRASH_DATE),
    total_injuries = INJURIES_TOTAL
  ) %>%
  group_by(year) %>%
  summarize(total_injuries = sum(total_injuries, na.rm = TRUE))

# Create the injuries visualization
ggplot(injuries_by_year, aes(x = year, y = total_injuries)) +
  geom_col(fill = "#9370DB", width = 0.7) +
  geom_text(
    aes(label = total_injuries),
    vjust = -0.5,
    size = 3.5,
    family = "Arial"
  ) +
  labs(
    x = "Year",
    y = NULL,
    title = "Injuries From Crashes Along DuSable Lake Shore Drive in District 6",
    subtitle = "Total persons sustaining fatal, incapacitating, non-incapacitating, and possible injuries",
    caption = "Data: Chicago City Data Portal \"Traffic Crashes\" Retrieved Feb. 6, 2025
Includes crashes reported within 100ft of DuSable Lake Shore Drive"
  ) +
  scale_x_continuous(
    breaks = seq(min_year, max_year, by = 1)
  ) +
  scale_y_continuous(
    expand = expansion(mult = c(0, 0.15)),
    breaks = function(x) seq(0, ceiling(max(x)), by = 50)
  ) +
  theme_classic() +
  theme(
    plot.title = element_text(face = "bold", size = 12, margin = margin(b = 10)),
    plot.subtitle = element_text(size = 10, color = "#666666"),
    plot.background = element_rect(fill = "#FFFCF7"),
    panel.background = element_rect(fill = "#FFFCF7"),
    panel.grid.major.y = element_line(color = "#f0f0f0", linewidth = 0.3),
    text = element_text(family = "Arial"),
    axis.line = element_line(color = "#666666"),
    axis.ticks = element_line(color = "#666666"),
    axis.title = element_text(color = "#666666"),
    axis.text = element_text(color = "#666666"),
    plot.caption = element_text(
      color = "#666666",
      size = 8,
      margin = margin(t = 10),
      hjust = 0
    ),
    plot.margin = margin(t = 15, r = 20, b = 15, l = 15, unit = "pt")
  )

Incapacitating Injuries Over Time

# Prepare incapacitating injuries data
incap_injuries_by_year <- lakeshore_crashes %>%
  st_drop_geometry() %>%
  mutate(
    year = lubridate::year(CRASH_DATE),
    incap_injuries = INJURIES_INCAPACITATING
  ) %>%
  group_by(year) %>%
  summarize(total_incap_injuries = sum(incap_injuries, na.rm = TRUE))

# Create the incapacitating injuries visualization
ggplot(incap_injuries_by_year, aes(x = year, y = total_incap_injuries)) +
  geom_col(fill = "#9370DB", width = 0.7) +
  geom_text(
    aes(label = total_incap_injuries),
    vjust = -0.5,
    size = 3.5,
    family = "Arial"
  ) +
  labs(
    x = "Year",
    y = NULL,
    title = "Incapacitating Injuries Along DuSable Lake Shore Drive in District 6",
    subtitle = "Total persons sustaining incapacitating/serious injuries in a crash",
    caption = "Data: Chicago City Data Portal \"Traffic Crashes\" Retrieved Feb. 6, 2025
Includes crashes reported within 100ft of DuSable Lake Shore Drive"
  ) +
  scale_x_continuous(
    breaks = seq(min_year, max_year, by = 1)
  ) +
  scale_y_continuous(
    expand = expansion(mult = c(0, 0.15)),
    breaks = function(x) seq(0, ceiling(max(x)), by = 5)
  ) +
  theme_classic() +
  theme(
    plot.title = element_text(face = "bold", size = 12, margin = margin(b = 10)),
    plot.subtitle = element_text(size = 10, color = "#666666"),
    plot.background = element_rect(fill = "#FFFCF7"),
    panel.background = element_rect(fill = "#FFFCF7"),
    panel.grid.major.y = element_line(color = "#f0f0f0", linewidth = 0.3),
    text = element_text(family = "Arial"),
    axis.line = element_line(color = "#666666"),
    axis.ticks = element_line(color = "#666666"),
    axis.title = element_text(color = "#666666"),
    axis.text = element_text(color = "#666666"),
    plot.caption = element_text(
      color = "#666666",
      size = 8,
      margin = margin(t = 10),
      hjust = 0
    ),
    plot.margin = margin(t = 15, r = 20, b = 15, l = 15, unit = "pt")
  )

Interactive Map Visualization

District 6 Fatality Map

Create an interactive map showing fatality locations in District 6:

# Transform geometries to WGS84 for mapping
lakeshore_fatalities <- st_transform(lakeshore_fatalities, 4326)
lake_shore_streets <- st_transform(lake_shore_streets, 4326)
buffer <- st_transform(buffer, 4326)
district_6 <- st_transform(district_6, 4326)

# Create Lake Shore Drive district 6 fatality map
district6_fatality_map <- leaflet() %>%
  addProviderTiles(providers$CartoDB.Positron) %>%
  # Add district 6
  addPolygons(
    data = district_6,
    fillColor = "#9370DB",  # Medium purple
    fillOpacity = 0.3,
    color = "#4B0082",  # Indigo
    weight = 2,
    label = ~paste("Senate District", senate_district)
  ) %>%
  # Add the buffer area
  addPolygons(
    data = buffer,
    fillColor = "#CCCCCC",  # Light gray
    fillOpacity = 0.2,
    color = "#CCCCCC",  # Light gray
    weight = 2
  ) %>%
  # Add Lake Shore Drive
  addPolylines(
    data = lake_shore_streets,
    color = "#CCCCCC",
    weight = 3
  ) %>%
  # Add clustered fatality circles
  addCircleMarkers(
    data = lakeshore_fatalities,
    lng = ~st_coordinates(geometry)[, 1],
    lat = ~st_coordinates(geometry)[, 2],
    radius = 6,  # Base size for individual points
    color = "red",
    fillColor = "red",
    fillOpacity = 0.6,
    stroke = TRUE,
    weight = 1,
    popup = ~paste(
      "Date:", CRASH_DATE, "<br>",
      "Fatalities:", fatality_count, "<br>",
      "Type:", FIRST_CRASH_TYPE
    ),
    clusterOptions = markerClusterOptions(
      maxClusterRadius = 50 * 0.0003,  # Approximate conversion of 50ft to degrees
      spiderfyOnMaxZoom = TRUE,
      iconCreateFunction = JS("
        function(cluster) {
          var count = cluster.getChildCount();
          var size = Math.min(50, 20 + Math.sqrt(count) * 8);
          return new L.DivIcon({
            html: '<div style=\"background-color: rgba(255, 0, 0, 0.6); width: ' + size + 'px; height: ' + size + 'px; margin-left: -' + (size/2) + 'px; margin-top: -' + (size/2) + 'px; border-radius: 50%; display: flex; align-items: center; justify-content: center; color: white; font-weight: bold; border: 2px solid rgba(255, 0, 0, 0.8);\">' + count + '</div>',
            className: 'marker-cluster-custom',
            iconSize: L.point(size, size)
          });
        }
      ")
    )
  ) %>%
  setView(lng = -87.6298, lat = 41.8781, zoom = 12)

# Display the map
district6_fatality_map

Lake Shore Drive Crash Maps

Create interactive maps showing injuries, incapacitating injuries, and fatalities along Lake Shore Drive:

# Transform geometries to WGS84 for mapping
lakeshore_crashes_4326 <- st_transform(lakeshore_crashes, 4326)
lakefront_districts_4326 <- st_transform(lakefront_districts, 4326)

# Create injuries heatmap
injuries_map <- leaflet() %>%
  addProviderTiles(providers$CartoDB.Positron) %>%
  # Add lakefront districts outlines
  addPolygons(
    data = lakefront_districts_4326,
    fillOpacity = 0,
    color = "#4B0082",  # Indigo
    weight = 2,
    label = ~paste("Senate District", senate_district)
  ) %>%
  # Add Lake Shore Drive
  addPolylines(
    data = lake_shore_streets,
    color = "#666666",
    weight = 3
  ) %>%
  # Add injuries heatmap
  addHeatmap(
    data = lakeshore_crashes_4326,
    lng = ~st_coordinates(geometry)[, 1],
    lat = ~st_coordinates(geometry)[, 2],
    intensity = ~INJURIES_TOTAL,
    blur = 20,
    max = 0.05,
    radius = 15
  ) %>%
  setView(lng = -87.6298, lat = 41.8781, zoom = 12)

# Create incapacitating injuries heatmap
incap_injuries_map <- leaflet() %>%
  addProviderTiles(providers$CartoDB.Positron) %>%
  # Add lakefront districts outlines
  addPolygons(
    data = lakefront_districts_4326,
    fillOpacity = 0,
    color = "#4B0082",  # Indigo
    weight = 2,
    label = ~paste("Senate District", senate_district)
  ) %>%
  # Add Lake Shore Drive
  addPolylines(
    data = lake_shore_streets,
    color = "#666666",
    weight = 3
  ) %>%
  # Add incapacitating injuries heatmap
  addHeatmap(
    data = lakeshore_crashes_4326,
    lng = ~st_coordinates(geometry)[, 1],
    lat = ~st_coordinates(geometry)[, 2],
    intensity = ~INJURIES_INCAPACITATING,
    blur = 20,
    max = 0.05,
    radius = 15
  ) %>%
  setView(lng = -87.6298, lat = 41.8781, zoom = 12)

# Create fatalities map with clusters
fatalities_map <- leaflet() %>%
  addProviderTiles(providers$CartoDB.Positron) %>%
  # Add lakefront districts outlines
  addPolygons(
    data = lakefront_districts_4326,
    fillOpacity = 0,
    color = "#4B0082",  # Indigo
    weight = 2,
    label = ~paste("Senate District", senate_district)
  ) %>%
  # Add Lake Shore Drive
  addPolylines(
    data = lake_shore_streets,
    color = "#666666",
    weight = 3
  ) %>%
  # Add clustered fatality circles
  addCircleMarkers(
    data = lakeshore_fatalities,
    lng = ~st_coordinates(geometry)[, 1],
    lat = ~st_coordinates(geometry)[, 2],
    radius = 6,
    color = "red",
    fillColor = "red",
    fillOpacity = 0.6,
    stroke = TRUE,
    weight = 1,
    popup = ~paste(
      "Date:", CRASH_DATE, "<br>",
      "Fatalities:", fatality_count, "<br>",
      "Type:", FIRST_CRASH_TYPE
    ),
    clusterOptions = markerClusterOptions(
      maxClusterRadius = 50 * 0.0003,
      spiderfyOnMaxZoom = TRUE,
      iconCreateFunction = JS("
        function(cluster) {
          var count = cluster.getChildCount();
          var size = Math.min(50, 20 + Math.sqrt(count) * 8);
          return new L.DivIcon({
            html: '<div style=\"background-color: rgba(255, 0, 0, 0.6); width: ' + size + 'px; height: ' + size + 'px; margin-left: -' + (size/2) + 'px; margin-top: -' + (size/2) + 'px; border-radius: 50%; display: flex; align-items: center; justify-content: center; color: white; font-weight: bold; border: 2px solid rgba(255, 0, 0, 0.8);\">' + count + '</div>',
            className: 'marker-cluster-custom',
            iconSize: L.point(size, size)
          });
        }
      ")
    )
  ) %>%
  setView(lng = -87.6298, lat = 41.8781, zoom = 12)

# Display the maps
injuries_map
incap_injuries_map
fatalities_map

Summary Statistics

# Calculate Lake Shore Drive statistics
lakeshore_total_crashes <- nrow(lakeshore_crashes)
lakeshore_total_injuries <- sum(lakeshore_crashes$INJURIES_TOTAL, na.rm = TRUE)
lakeshore_total_incap_injuries <- sum(lakeshore_crashes$INJURIES_INCAPACITATING, na.rm = TRUE)
lakeshore_total_fatalities <- sum(!is.na(lakeshore_crashes$FATALITY_VICTIM) & lakeshore_crashes$FATALITY_VICTIM != "", na.rm = TRUE)

# Calculate total district statistics
district_total_crashes <- nrow(crashes_with_district)
district_total_injuries <- sum(crashes_with_district$INJURIES_TOTAL, na.rm = TRUE)
district_total_incap_injuries <- sum(crashes_with_district$INJURIES_INCAPACITATING, na.rm = TRUE)
district_total_fatalities <- sum(!is.na(crashes_with_district$FATALITY_VICTIM) & crashes_with_district$FATALITY_VICTIM != "", na.rm = TRUE)

# Calculate percentages
pct_crashes <- (lakeshore_total_crashes / district_total_crashes) * 100
pct_injuries <- (lakeshore_total_injuries / district_total_injuries) * 100
pct_incap_injuries <- (lakeshore_total_incap_injuries / district_total_incap_injuries) * 100
pct_fatalities <- (lakeshore_total_fatalities / district_total_fatalities) * 100

cat("Lake Shore Drive District 6 Analysis Summary\n")
## Lake Shore Drive District 6 Analysis Summary
cat("----------------------------------------\n")
## ----------------------------------------
cat("Lake Shore Drive Statistics:\n")
## Lake Shore Drive Statistics:
cat("Total Crashes:", lakeshore_total_crashes, "\n")
## Total Crashes: 6248
cat("Total Injuries:", lakeshore_total_injuries, "\n")
## Total Injuries: 1415
cat("Total Incapacitating Injuries:", lakeshore_total_incap_injuries, "\n")
## Total Incapacitating Injuries: 166
cat("Total Fatalities:", lakeshore_total_fatalities, "\n\n")
## Total Fatalities: 16
cat("Percentage of District 6 Totals:\n")
## Percentage of District 6 Totals:
cat(sprintf("Crashes: %.1f%%\n", pct_crashes))
## Crashes: 19.2%
cat(sprintf("Injuries: %.1f%%\n", pct_injuries))
## Injuries: 26.0%
cat(sprintf("Incapacitating Injuries: %.1f%%\n", pct_incap_injuries))
## Incapacitating Injuries: 28.4%
cat(sprintf("Fatalities: %.1f%%\n", pct_fatalities))
## Fatalities: 61.5%